# Loading required package
library(margins)
# Loading the data
load("analysis_1_data.RData")

#######################################
# IN-TEXT DESCRIPTION OF MEAN SUPPORT #
#######################################
# Overall
round(mean(nominees$support_delay, na.rm=T),3)*100

# Co-partisans of president (including leaners)
round(mean(nominees$support_delay[which(nominees$party_aligned_with_president == 1)], na.rm=T),3)*100

# Out-partisans of president (including leaners)
round(mean(nominees$support_delay[which(nominees$outpartisan_of_president == 1)], na.rm=T),3)*100

# Independents (true independents)
round(mean(nominees$support_delay[which(nominees$independent == 1)], na.rm=T),3)*100

############
# FIGURE 1 #
############
# Removing independents 
nominees.partisans <- nominees[which(nominees$independent == 0),]

# Baseline model; results discussed in text (Model 1 in Appendix table C.3): both nominations, no interactions
model.1 <- lm(support_delay ~ party_aligned_with_president
              + sc_approve + knowledge_of_vacancy + white + age_numeric + male 
              + education + scalia_vacancy, weights=weights, data=nominees.partisans)
summary(model.1)

# Model used for building the figure (Model 4 in Appendix table C.3): Pooling both nominations, triple interaction
model.4.triple.interaction <- lm(support_delay ~ party_aligned_with_president*sc_approve*knowledge_of_vacancy 
                                 + party_aligned_with_president + knowledge_of_vacancy + sc_approve + white + age_numeric 
                                 + male + education + scalia_vacancy, weights=weights, data=nominees.partisans)

### Calculating marginal effects ###
# Marginal effect of SC approve by knowledge and partisan alignment with president
margins.of.model.4.approval.both <- margins(model.4.triple.interaction, variables='sc_approve', at=list(knowledge_of_vacancy=c(0,1), party_aligned_with_president=c(0,1)))
summary(margins.of.model.4.approval.both)

# Marginal effect of knowledge by SC approve and partisan alignment with president
margins.of.model.4.knowledge.both <- margins(model.4.triple.interaction, variables='knowledge_of_vacancy', at=list(sc_approve=c(0,1), party_aligned_with_president=c(0,1)))
summary(margins.of.model.4.knowledge.both)

# Marginal effect of SC approval, averaging across knowledge
margins.of.model.4.approval <- margins(model.4.triple.interaction, variables='sc_approve', at=list(party_aligned_with_president=c(0,1)))
summary(margins.of.model.4.approval)

# Marginal effect of knowledge, averaging across approval
margins.of.model.4.knowledge <- margins(model.4.triple.interaction, variables='knowledge_of_vacancy', at=list(party_aligned_with_president=c(0,1)))
summary(margins.of.model.4.knowledge)

### Left panel -- results for approval ###
# Averaged across knowledge
# Effect, outpartisans
outpartisans <- c(summary(margins.of.model.4.approval)[1,"AME"])
# Confidence interval
outpartisans.lower <- c(summary(margins.of.model.4.approval)[1,"lower"])
outpartisans.upper <- c(summary(margins.of.model.4.approval)[1,"upper"])

# Effect, copartisans
copartisans <- c(summary(margins.of.model.4.approval)[2,"AME"])
# Confidence interval
copartisans.lower <- c(summary(margins.of.model.4.approval)[2,"lower"])
copartisans.upper <- c(summary(margins.of.model.4.approval)[2,"upper"])

# Separate for low and high knowledge
# Effect
outpartisans.approval <- c(summary(margins.of.model.4.approval.both)[1,"AME"],summary(margins.of.model.4.approval.both)[3,"AME"])
copartisans.approval <- c(summary(margins.of.model.4.approval.both)[2,"AME"],summary(margins.of.model.4.approval.both)[4,"AME"])

# Confidence interval
outpartisans.approval.lower <- c(summary(margins.of.model.4.approval.both)[1,"lower"],summary(margins.of.model.4.approval.both)[3,"lower"])
copartisans.approval.lower <- c(summary(margins.of.model.4.approval.both)[2,"lower"],summary(margins.of.model.4.approval.both)[4,"lower"])
outpartisans.approval.upper <- c(summary(margins.of.model.4.approval.both)[1,"upper"],summary(margins.of.model.4.approval.both)[3,"upper"])
copartisans.approval.upper <- c(summary(margins.of.model.4.approval.both)[2,"upper"],summary(margins.of.model.4.approval.both)[4,"upper"])

# Plotting coordinates
group <- c(0.8,3.2)
group2 <- c(0.1,1.5,2.5,3.9)
group3 <- c(0.1,0.8,1.5,2.5,3.2,3.9)

par(mar=c(4.5,5,1,1), mgp=c(2,.25,0))
plot(group,c(outpartisans,copartisans),xlim=c(-0.25,4.25),ylim=c(-0.3,0.2), tck=0, cex.main=1, main="", cex.lab=1.3, cex.axis=1.2, 
     pch=c(16,17), cex=2.2, ylab="Marginal Effect of SC Approval on Support for Delay", xaxt="n", yaxt="n", xlab="", type="p")
axis(2,at=seq(-0.3,0.2,by=.1), cex.axis=1, 
     labels=c("-0.3","-0.2","-0.1","0","0.1","0.2"), tck=0.01, las=1)
axis(1, at=group3,cex.axis=0.8,labels=c("No Knowledge","Average","Knowledge","No Knowledge","Average","Knowledge"), tck=-0.01)
axis(1, line=1.5, at=group,cex.axis=1.2, labels=c("Outpartisans","Copartisans"), lty="blank")
points(group2, c(outpartisans.approval,copartisans.approval), pch=c(16,16,17,17), col="grey30", cex=1.6)
segments(group, c(outpartisans.lower,copartisans.lower), group, c(outpartisans.upper,copartisans.upper), lwd=2.75)
segments(group2, c(outpartisans.approval.lower,copartisans.approval.lower), 
         group2, c(outpartisans.approval.upper,copartisans.approval.upper), col='grey30', lwd=2.25, lty=2)
abline(h=0)
# Plotting text
# Outpartisans, averaged
text(0.8, -0.10, "-0.17 (p<0.01)", cex=1.35)
# Copartisans, averaged
text(3.2, -0.1, "-0.02 (p<0.49)", cex=1.35)
# Outpartisans, by knowledge
text(0.1, -0.165, "-0.06 (p<0.21)", cex=1.05, col="grey30")
text(1.5, -0.295, "-0.22 (p<0.01)", cex=1.05, col="grey30")
# Copartisans, by knowledge
text(2.5, 0.13, "-0.01 (p<0.90)", cex=1.05, col="grey30")
text(3.9, 0.06, "-0.02 (p<0.45)", cex=1.05, col="grey30")

########################################
# Right panel -- results for knowledge #
########################################
# Averaged across approval
# Effect, outpartisans
outpartisans <- c(summary(margins.of.model.4.knowledge)[1,"AME"])
# Confidence interval
outpartisans.lower <- c(summary(margins.of.model.4.knowledge)[1,"lower"])
outpartisans.upper <- c(summary(margins.of.model.4.knowledge)[1,"upper"])

# Effect, copartisans
copartisans <- c(summary(margins.of.model.4.knowledge)[2,"AME"])
# Confidence interval
copartisans.lower <- c(summary(margins.of.model.4.knowledge)[2,"lower"])
copartisans.upper <- c(summary(margins.of.model.4.knowledge)[2,"upper"])

# Separate for low and high SC approval
# Effect
outpartisans.knowledge <- c(summary(margins.of.model.4.knowledge.both)[1,"AME"],summary(margins.of.model.4.knowledge.both)[3,"AME"])
copartisans.knowledge <- c(summary(margins.of.model.4.knowledge.both)[2,"AME"],summary(margins.of.model.4.knowledge.both)[4,"AME"])

# Confidence interval
outpartisans.knowledge.lower <- c(summary(margins.of.model.4.knowledge.both)[1,"lower"],summary(margins.of.model.4.knowledge.both)[3,"lower"])
copartisans.knowledge.lower <- c(summary(margins.of.model.4.knowledge.both)[2,"lower"],summary(margins.of.model.4.knowledge.both)[4,"lower"])
outpartisans.knowledge.upper <- c(summary(margins.of.model.4.knowledge.both)[1,"upper"],summary(margins.of.model.4.knowledge.both)[3,"upper"])
copartisans.knowledge.upper <- c(summary(margins.of.model.4.knowledge.both)[2,"upper"],summary(margins.of.model.4.knowledge.both)[4,"upper"])

# Plotting coordinates
group <- c(0.8,3.2)
group2 <- c(0.1,1.5,2.5,3.9)
group3 <- c(0.1,0.8,1.5,2.5,3.2,3.9)

par(mar=c(4.5,5,1,1), mgp=c(2,.25,0))
plot(group,c(outpartisans,copartisans),xlim=c(-0.25,4.25),ylim=c(-0.3,0.2), tck=0, cex.main=1, main="", cex.lab=1.3, cex.axis=1.2, 
     pch=c(16,17), cex=2.4, ylab="Marginal Effect of Vacancy Knowledge on Support for Delay", xaxt="n", yaxt="n", xlab="", type="p")
axis(2,at=seq(-0.3,0.2,by=.1), cex.axis=1, 
     labels=c("-0.3","-0.2","-0.1","0","0.1","0.2"), tck=0.01, las=1)
axis(1, at=group3,cex.axis=0.8,labels=c("Don't Approve","Average","Approve","Don't Approve","Average","Approve"), tck=-0.01)
axis(1, line=1.5, at=group,cex.axis=1.2, labels=c("Outpartisans","Copartisans"), lty="blank")
points(group2, c(outpartisans.knowledge,copartisans.knowledge), pch=c(16,16,17,17), col="grey30", cex=1.6)
segments(group, c(outpartisans.lower,copartisans.lower), group, c(outpartisans.upper,copartisans.upper), lwd=2.75)
segments(group2, c(outpartisans.knowledge.lower,copartisans.knowledge.lower), 
         group2, c(outpartisans.knowledge.upper,copartisans.knowledge.upper), col='grey30', lwd=2.25, lty=2)
abline(h=0)
# Plotting text
# Outpartisans, averaged
text(0.8, -0.05, "0.02 (p<0.41)", cex=1.35)
# Copartisans, averaged
text(3.2, -0.22, "-0.13 (p<0.01)", cex=1.35)
# Outpartisans, by approval
text(0.1, .03, "0.13 (p<0.01)", cex=1.05, col="grey30")
text(1.5, 0.05, "-0.03 (p<0.36)", cex=1.05, col="grey30")
# Copartisans, by approval
text(2.5, -0.26, "-0.12 (p<0.03)", cex=1.05, col="grey30")
text(3.9, -0.05, "-0.14 (p<0.01)", cex=1.05, col="grey30")
